16.7 Энтропия и семейство экспоненциальных распределений

Энтропия

Информативность наблюдений

Представим, что вы пронаблюдали значение xx некоторой случайной величины ξ\xi. Как измерить количество информации h(x)h(x), которое вы при этом получили? Следующие соображения кажутся в этом плане вполне естественными:

  • чем выше вероятность P(ξ=x)\mathbb P(\xi = x), тем более ожидаемо появление значения xx и, соответственно, менее информативно;
  • и наоборот, наблюдение маловероятного значения xx обычно даёт обильную пищу для размышлений и повышает h(x)h(x);
  • при наблюдении двух независимых реализаций xx и yy случайной величины ξ\xi логично складывать полученную информацию: h((x,y))=h(x)+h(y)h((x, y)) = h(x) + h(y).

Указанные соображения наводят на мысль, что информацию следует считать убывающей функцией от вероятности: h(x)=g(P(ξ=x))h(x) = g\big(\mathbb P(\xi = x)\big). Кроме того, функция gg должна превращать произведение в сумму, поскольку для независимых случайных величин ξ\xi и η\eta
равенство h((x,y))=h(x)+h(y)h((x, y)) = h(x) + h(y) влечёт

g(P(ξ=x,η=y))=g(P(ξ=x)P(η=y))=g(P(ξ=x))+g(P(η=y)). g\big(\mathbb P(\xi = x, \eta=y)\big) = g\big(\mathbb P(\xi = x) \mathbb P(\eta = y)\big) = g\big(\mathbb P(\xi = x)\big) + g\big(\mathbb P(\eta = y)\big).

На самом деле выбор тут небогат. Единственная непрерывная функция, обладающая такими свойствами, — это логарифм: g(p)=logpg(p) = -\log p. Основание логарифма может быть любым числом больше единицы. Поскольку информацию измеряют в битах и байтах, в теории информации обычно предпочитают двоичные логарифмы. Однако для вычислений удобнее использовать натуральный логарифм, и по умолчанию мы будем подразумевать именно его (кстати, соответствующую единицу информации называют «нат»).

Энтропия Шеннона

Среднее количество информации, которое несёт в себе значение дискретной случайной ξ\xi с распределением вероятностей pk=P(ξ=k)p_k= \mathbb P(\xi = k), вычисляется по формуле

Hξ=E(g(p(ξ)))=E(logp(ξ))=kpklogpk. \mathbb H\xi = \mathbb E \big(g(p(\xi))\big) = -\mathbb E \big(\log p(\xi)\big)=-\sum\limits_k p_k \log p_k.

Это так называемая энтропия (Шэннона).

Небольшое математическое замечание

При вычислении энтропии регулярно можно встретить выражение 0log00\cdot \log 0 с не вполне очевидным значением. Поскольку limt0+tlogt=0\lim \limits_{t\to 0+} t\log t = 0, по определению полагаем 0log0=0.0\cdot \log 0 = 0.

Пример. Рассмотрим схему Бернулли с вероятностью «успеха» pp. Энтропия её результата равна

Hξ=(1p)log(1p)plogp,ξBern(p).\mathbb H\xi = -(1 - p)\log(1 - p) - p\log p, \quad \xi \sim \mathrm{Bern}(p).

Давайте посмотрим на график этой функции:

entropy

Минимальное значение (нулевое) энтропия принимает при p=0p = 0 или p=1p=1. Исход такого вырожденного эксперимента заранее известен, и чтобы сообщить кому-то о его результате, достаточно 00 бит информации. Иначе говоря, можно вообще ничего не передавать, и так всё предельно ясно.

Максимальное значение энтропии достигается в точке 12\frac12, что вполне соответствует тому, что при p=12p=\frac12 предсказать исход эксперимента сложнее всего.

Упражнение. Найдите энтропию геометрического распределения с вероятностью «успеха» pp: ξGeom(p)\xi \sim \mathrm{Geom}(p),

P(ξ=k)=p(1p)k1,kN,0<p1.\mathbb P(\xi = k) = p(1-p)^{k-1}, \quad k \in \mathbb N, \quad 0 < p \leqslant 1.

Решение (не открывайте сразу, попробуйте сначала решить самостоятельно)

По определению энтропии имеем

Hξ=k=1p(1p)k1(logp+(k1)log(1p))=\mathbb H\xi = -\sum\limits_{k=1}^\infty p(1-p)^{k-1} (\log p + (k-1)\log(1-p)) =

=logpplog(1p)k=1k(1p)k.= -\log p -p\log(1-p)\sum\limits_{k=1}^\infty k(1-p)^k.

Дифференцированием геометрической прогрессии находим

k=1kxk1=k=1(xk)=ddxk=0xk=(11x)=1(1x)2;\sum\limits_{k=1}^\infty kx^{k-1} = \sum\limits_{k=1}^\infty (x^k)' = \frac d{dx}\sum\limits_{k=0}^\infty x^k = \Big(\frac 1{1-x}\Big)' = \frac 1{(1-x)^2};

подставляя сюда x=1px = 1-p, окончательно получаем, что

Hξ=logpp(1p)log(1p)1p2=logp(1p)log(1p)p. \mathbb H\xi = -\log p -p(1-p)\log(1-p)\cdot \frac 1{p^2} = -\log p - \frac{(1-p)\log(1-p)}p.

Вспомним, что случайная величина ξGeom(p)\xi\sim\mathrm{Geom}(p) равна количеству независимых испытаний с вероятностью «успеха» pp до появления первого «успеха». При p=1p=1 энтропия Hξ\mathbb H\xi минимальна и равна нулю, ведь в этом случае геометрическое распределение становится вырожденным: «успех» наступает сразу же с вероятностью 11. А вот при p0+p\to 0+ серия неудачных испытаний может быть сколь угодно длинной; распределение становится всё более неопределённым и «размазанным» по своему носителю, и его энтропия
стремится к ++\infty. График энтропии как функции от pp выглядит так:

entropy

Следующие свойства энтропии дискретной случайной величины ξ\xi вытекают прямо из определения:

  • неотрицательность: Hξ0\mathbb H\xi \geqslant 0;
  • Hξ=0    P(ξ=a)=1\mathbb H\xi = 0 \iff \mathbb P(\xi = a) = 1 при некотором aRa\in\mathbb R (нулевую энтропию имеют вырожденные распределения и только они);
  • Hξlogn\mathbb H\xi \leqslant \log n, если случайная величина имеет конечный носитель мощности nn.

Последнее свойство выводится из неравенства Йенсена. Применяя его к выпуклой вверх логарифмической функции, с учётом нормировки условия k=1npk=1\sum\limits_{k=1}^n p_k = 1 получаем

k=1npklogpk=k=1npklog1pklog(k=1npk1pk)=logn. -\sum\limits_{k=1}^n p_k\log p_k = \sum\limits_{k=1}^n p_k\log \frac 1p_k \leqslant \log\bigg(\sum\limits_{k=1}^n p_k\cdot \frac 1{p_k}\bigg) = \log n.

Вопрос на подумать. Итак, всякое распределение с носителем {1,2,,n}\{1, 2, \ldots, n\} имеет энтропию не больше logn\log n. А у какого распределения она в точности равна logn\log n?

Ответ (не открывайте сразу, сначала подумайте самостоятельно)

Такое произойдёт, если все вероятности pk=1np_k = \frac 1n, 1kn1\leqslant k \leqslant n, т.е. распределение равномерное.

Дифференциальная энтропия

Чтобы вычислить энтропию непрерывной случайной величины ξ\xi, надо, как водится, сумму заменить на интеграл, и получится формула дифференциальной энтропии:

Hξ=pξ(x)logpξ(x)dx. \mathbb H\xi = -\int p_{\xi}(x) \log p_{\xi}(x)\, dx.

Замечание. В дальнейшем мы будем использовать одинаковый термин энтропия как для дискретных, так и для непрерывных случайных величин, для краткости опуская слово дифференциальная в последнем случае. Кроме того, энтропию распределения pp, заданного через pmf или pdf, будем обозначать H[p]\mathbb H[p]. Такое обозначение позволяет избежать привязки к случайной величине там, где это излишне. Если ξp(x)\xi \sim p(x), то обозначения Hξ\mathbb H\xi и H[p]\mathbb H[p] эквивалентны. Также отметим, что энтропию можно записать в виде математического ожидания:

H[p]=Eξp(x)log1p(ξ). \mathbb H[p] = \mathbb E_{\xi \sim p(x)} \log\frac 1{p(\xi)}.

Пример. Найдём энтропию нормального распределения N(μ,σ2)\mathcal N(\mu, \sigma^2). Его плотность равна p(x)=12πσe(xμ)22σ2p(x) = \frac 1{\sqrt{2\pi} \sigma} e^{-\frac{(x-\mu)^2}{2\sigma^2}}, следовательно,

H[p]=+p(x)(12ln(2πσ2)+(xμ)22σ2)dx=\mathbb H[p] = \int\limits_{-\infty}^{+\infty}p(x)\Big(\frac 12\ln\big(2\pi\sigma^2\big) + \frac{(x-\mu)^2}{2\sigma^2}\Big)dx =

=12ln(2πσ2)+12πσ+(xμ)22σ2e(xμ)22σ2dx.= \frac 12\ln\big(2\pi\sigma^2\big) + \frac 1{\sqrt{2\pi} \sigma}\int\limits_{-\infty}^{+\infty}\frac{(x-\mu)^2}{2\sigma^2}e^{-\frac{(x-\mu)^2}{2\sigma^2}}\,dx.

Делая в последнем интеграле замену t=(xμ)22σ2t = -\frac{(x-\mu)^2}{2\sigma^2}, получаем, что

H[p]=12ln(2πσ2)+12π0+2tetdt=12ln(2πσ2)+Γ(32)π.\mathbb H[p] = \frac 12\ln\big(2\pi\sigma^2\big) + \frac 1{\sqrt{2\pi}}\int\limits_{0}^{+\infty}\sqrt{2t}e^{-t}\,dt =\frac 12\ln\big(2\pi\sigma^2\big) + \frac{\Gamma\big(\frac 32\big)}{\sqrt \pi}.

По свойству гамма-функции Γ(32)=12Γ(12)=π2\Gamma\big(\frac 32\big) = \frac 12 \Gamma\big(\frac 12\big) = \frac{\sqrt \pi}2. Таким образом,

H[p]=12ln(2πσ2)+12. \mathbb H[p] = \frac 12\ln\big(2\pi\sigma^2\big) + \frac 12.

Как видно, энтропия гауссовского распределения N(μ,σ2)\mathcal N(\mu, \sigma^2) не ограничена ни сверху, ни снизу:

limσ+12ln(2πσ2)=+,limσ0+12ln(2πσ2)=. \lim\limits_{\sigma\to +\infty} \frac 12\ln\big(2\pi\sigma^2\big) = +\infty, \quad \lim\limits_{\sigma\to 0+} \frac 12\ln\big(2\pi\sigma^2\big) = -\infty.

И да, в отличие от энтропии дискретного распределения дифференциальная энтропия может быть отрицательной. Это связано с тем, что плотность может принимать значения больше единицы, и поэтому математическое ожидание её логарифма с обратным знаком может оказаться меньше нуля. В частности, с нормальным распределением так происходит, если σ2<12πe\sigma^2 < \frac 1{2\pi e}.

Упражнение. Найдите энтропию показательного распределения Exp(λ)\mathrm{Exp}(\lambda).

Решение (не открывайте сразу, попробуйте сначала решить самостоятельно)

Плотность случайной величины ξExp(λ)\xi \sim \mathrm{Exp}(\lambda) равна pξ(x)=λeλxp_{\xi}(x) = \lambda e^{-\lambda x}, x0x \geqslant 0. Таким образом,

Hξ=H[pξ]=0+pξ(x)(λxlogλ)dx=0+λ2xeλxdxlogλ. \mathbb H\xi = \mathbb H[p_\xi] = \int\limits_{0}^{+\infty} p_{\xi}(x)(\lambda x - \log \lambda)\, dx = \int\limits_{0}^{+\infty} \lambda^2 xe^{-\lambda x}\,dx - \log\lambda.

С помощью замены t=λxt = \lambda x находим, что

Hξ=0+tetdtlogλ=Γ(2)logλ=1logλ. \mathbb H\xi = \int\limits_{0}^{+\infty} te^{-t}\,dt - \log\lambda = \Gamma(2) - \log\lambda = 1 - \log\lambda.

KL-дивергенция

В задачах машинного обучения истинное распределение p(x)p(x), из которого приходят наблюдения, обычно неизвестно, и его пытаются приблизить распределением q(x)q(x) из некоторого класса модельных распределений. Дивергенция Кульбака-Лейблера (KL-дивергенция, относительная энтропия) позволяет оценить расстояние между распределениями pp и qq:

KL(pq)=kpklogpkqk\mathbb{KL}(p\vert\vert q) = \sum\limits_k p_k\log\frac{p_k}{q_k}

в дискретном случае и

KL(pq)=p(x)logp(x)q(x)dx\mathbb{KL}(p\vert\vert q) = \int p(x)\log \frac{p(x)}{q(x)}dx

в непрерывном. KL-дивергенцию можно представить в виде разности:

KL(pq)=p(x)log1q(x)dxp(x)log1p(x)dx=\mathbb{KL}(p\vert\vert q) = \int p(x)\log \frac 1{q(x)}dx - \int p(x)\log\frac 1{p(x)}dx =

=Eξp(x)1logq(ξ)кросс-энтропияEξp(x)1logp(ξ)энтропия.=\underbrace{\mathbb E_{\xi \sim p(x)} \frac 1{\log q(\xi)}}_{\text{кросс-энтропия}} - \underbrace{\mathbb E_{\xi \sim p(x)} \frac1{\log p(\xi)}}_{\text{энтропия}}.

Здесь вычитаемое – это уже знакомая нам энтропия распределения p(x)p(x), которая показывает, сколько в среднем бит требуется, чтобы закодировать значение случайной величины ξp(x)\xi \sim p(x). Уменьшаемое носит название кросс-энтропии распределений p(x)p(x) и q(x)q(x).
Кросс-энтропию можно интерпретировать как среднее число бит для кодирования значения случайной величины ξp(x)\xi \sim p(x) алгоритмом, оптимизированным для кодирования случайной величины ηq(x)\eta \sim q(x). Иными словами, дивергенция Кульбака-Лейблера говорит о том, насколько увеличится средняя длина кодов для значений pp, если при настройке алгоритма кодирования вместо pp использовать qq. Подробнее об этом вы можете почитать, например, в данном посте.

Дивергенция Кульбака-Лейблера в некотором роде играет роль расстояния между распределениями. В частности, KL(pq)0\mathbb{KL}(p\vert\vert q)\geqslant 0, причём дивергенция равна нулю, только если распределения совпадают почти всюду (для дискретных и непрерывных распределений это означает, что они просто тождественны). Но при этом она не является симметричной: вообще говоря, KL(pq)KL(qp)\mathbb{KL}(p\vert\vert q)\ne \mathbb{KL}(q\vert\vert p).

Упражнение. Пользуясь неравенством ln(1+t)t\ln (1+t) \leqslant t, t>1t > -1, докажите неотрицательность KL-дивергенции.

Попробуйте доказать самостоятельно, прежде чем смотреть решение

Указанное неравенство можно переписать в виде lntt1-\ln t \geqslant t - 1, t>0t > 0, поэтому

KL(pq)=p(x)lnq(x)p(x)dx \mathbb{KL}(p\vert\vert q) = -\int p(x)\ln{\frac{q(x)}{p(x)}}dx \geqslant

p(x)(q(x)p(x)1)dx=q(x)dxp(x)dx=0. \geqslant \int p(x)\Big(\frac{q(x)}{p(x)} - 1\Big)dx = \int q(x)dx - \int p(x)dx = 0.

Пример. С помощью KL-дивергенции измерим расстояние между двумя гауссианами p(x)=N(xμ1,σ12)p(x) = \mathcal N(x \vert \mu_1, \sigma_1^2) и q(x)=N(xμ2,σ22)q(x) = \mathcal N(x \vert \mu_2, \sigma_2^2).

Подставляя явные выражения для плотностей

p(x)=12πσ1e(xμ1)22σ12 и q(x)=12πσ2e(xμ2)22σ22,p(x) = \frac 1{\sqrt{2\pi}\sigma_1} e^{-\frac{(x-\mu_1)^2}{2\sigma_1^2}}\text{ и } q(x) = \frac 1{\sqrt{2\pi}\sigma_2} e^{-\frac{(x-\mu_2)^2}{2\sigma_2^2}},

находим

lnp(x)q(x)=lnσ2σ1+(xμ2)22σ22(xμ1)22σ12=\ln \frac{p(x)}{q(x)} = \ln \frac{\sigma_2}{\sigma_1} + \frac{(x-\mu_2)^2}{2\sigma_2^2} -\frac{(x-\mu_1)^2}{2\sigma_1^2} =

=lnσ2σ1+(12σ2212σ12)(xμ1)2+μ1μ2σ22(xμ1)+(μ1μ2)22σ22.=\ln \frac{\sigma_2}{\sigma_1} + \Big(\frac 1{2\sigma_2^2} - \frac 1{2\sigma_1^2}\Big)(x-\mu_1)^2 + \frac{\mu_1 - \mu_2}{\sigma_2^2}(x-\mu_1) + \frac{(\mu_1 - \mu_2)^2}{2\sigma_2^2}.

Из свойств нормального распределения вытекает, что

+p(x)(xμ1)dx=0,+p(x)(xμ1)2dx=σ12. \int\limits_{-\infty}^{+\infty}p(x) (x-\mu_1)\,dx = 0, \quad \int\limits_{-\infty}^{+\infty}p(x) (x-\mu_1)^2\,dx = \sigma_1^2.

Таким образом,

KL(pq)=Eξp(x)lnp(ξ)q(ξ)=lnσ2σ1+σ12+(μ1μ2)22σ2212.\mathbb{KL}(p\vert\vert q) = \mathbb E_{\xi \sim p(x)} \ln \frac{p(\xi)}{q(\xi)} = \ln \frac{\sigma_2}{\sigma_1} + \frac{\sigma_1^2 + (\mu_1 - \mu_2)^2}{2\sigma_2^2} - \frac 12.

Как и должно быть, полученное выражение равно нулю, если гауссианы совпадают. При равных дисперсиях σ12=σ22=σ2\sigma_1^2 = \sigma_2^2 = \sigma^2 получаем, что KL(pq)=(μ1μ2)22σ2\mathbb{KL}(p\vert\vert q) = \frac{(\mu_1 - \mu_2)^2}{2\sigma^2}. Это выражение остаётся прежним, если поменять местами μ1\mu_1 и μ2\mu_2, поэтому в этом случае KL(pq)=KL(qp)\mathbb{KL}(p\vert\vert q) = \mathbb{KL}(q\vert\vert p). Если же σ1σ2\sigma_1 \ne \sigma_2, то выражение для KL(qp)\mathbb{KL}(q\vert\vert p) явно отличается от KL(pq)\mathbb{KL}(p\vert\vert q), что лишний раз показывает несимметричность KL-дивергенции.

Упражнение. Найдите дивергенцию Кульбака-Лейблера двух показательных распределений p(x)=Exp(xλ)p(x) = \mathrm{Exp}(x \vert \lambda) и q(x)=Exp(xμ)q(x) = \mathrm{Exp}(x \vert \mu).

Попробуйте решить самостоятельно, прежде чем смотреть решение

Имеем p(x)=λeλxp(x) = \lambda e^{-\lambda x}, q(x)=μeμxq(x) = \mu e^{-\mu x}, x0x\geqslant 0, откуда

lnp(x)q(x)=lnλμ+(μλ)x.\ln \frac{p(x)}{q(x)} = \ln \frac \lambda\mu + (\mu-\lambda)x.

Следовательно,

KL(pq)=lnλμ+(μλ)0+λxeλxdx= \mathbb{KL}(p\vert\vert q) = \ln \frac \lambda\mu + (\mu-\lambda)\int\limits_0^{+\infty} \lambda x e^{-\lambda x}\,dx =

=lnλμ+μλλ0+tetdt=lnλμ+μλ1. =\ln\frac \lambda\mu + \frac{\mu-\lambda}{\lambda}\int\limits_0^{+\infty} te^{-t}dt = \ln \frac \lambda\mu + \frac \mu \lambda - 1.

И здесь KL-дивергенция равна нулю при λ=μ\lambda = \mu.

Кросс-энтропия

При определении KL-дивергенции мы уже встречались с кросс-энтропией

H[p,q]=Eξp(x)logq(ξ) \mathbb H[p, q] = -\mathbb E_{\xi \sim p(x)} \log q(\xi)

В зависимости от типа распределений кросс-энтропия вычисляется по формуле

H[p,q]=kpklogqk или H[p,q]=+p(x)logq(x)dx. \mathbb H[p, q] = -\sum\limits_k p_k \log q_k \text{ или } \mathbb H[p, q] = -\int\limits_{-\infty}^{+\infty} p(x) \log q(x)\, dx.

Поскольку

KL(pq)=H[p,q]H[p], \mathbb{KL}(p\vert\vert q) = \mathbb H[p, q] - \mathbb H[p],

задача минимизации KL-дивергенции между неизвестным распределением данных p(x)p(x) и модельным распределением q(x)q(x) эквивалентна задаче минимизации кросс-энтропии. Разница между ними равна энтропии распределения p(x)p(x), которая, очевидно, не зависит от q(x)q(x).

В машинном обучении кросс-энтропию часто используют в качестве функции потерь в задаче классификации на K>1K>1 классов. Истинное распределение на каждом обучающем объекте задаётся с помощью one hot encoding и является вырожденным:

y=(y1,,yK),yk{0,1},k=1Kyk=1.\boldsymbol y = (y_1, \ldots, y_K), \quad y_k \in \{0,1\}, \quad \sum\limits_{k=1}^K y_k = 1.

Классификатор обычно выдаёт вероятности принадлежности каждому из классов,

y^=(y^1,,y^K),y^k=P(класс k).\boldsymbol{\widehat y} = (\widehat y_1, \ldots, \widehat y_K), \quad \widehat y_k = \mathbb P(\text{класс }k).

Функция потерь на одном объекте полагается равной кросс-энтропии между истинным и предсказанным распределениями:

L(y,y^)=k=1Kyklogy^k. \mathcal L(\boldsymbol y, \boldsymbol {\widehat y}) = -\sum\limits_{k=1}^K y_k \log \widehat y_k.

И это вполне логично: чем ближе модельное распределение к истинному, тем меньше наши потери. В идеале L(y,y^)=0\mathcal L(\boldsymbol y, \boldsymbol {\widehat y}) = 0, если y=y^\boldsymbol y = \boldsymbol {\widehat y}.

Чтобы вычислить функцию потерь по обучающей выборке из NN объектов с метками y(i)\boldsymbol y^{(i)}, обычно берут усреднённную кросс-энтропию

1Ni=1NL(y(i),y^(i))=1Ni=1Nk=1Kyk(i)logy^k(i). \frac 1N \sum\limits_{i=1}^N\mathcal L(\boldsymbol y^{(i)}, \boldsymbol{\widehat y}^{(i)}) =-\frac 1N\sum\limits_{i=1}^N\sum\limits_{k=1}^K y^{(i)}_k \log \widehat y^{(i)}_k.

Принцип максимальной энтропии

В параграфе про оценки параметров были описаны различные свойства параметрических оценок и методы их получения, например, метод моментов или метод максимального правдоподобия. В принципе, если мы уже выбрали для наших данных X1,,XnX_1, \ldots, X_n некоторое параметрическое семейство Fθ(x)F_\theta(x), моделирующее их распределение, восстановить его параметры чаще всего можно по выборочному среднему Xn=1nk=1nXk\overline{X}_n = \frac{1}{n}\sum\limits_{k = 1}^n X_k и/или выборочной дисперсии Sn=1nk=1n(XkXn)2\overline S_n = \frac{1}{n} \sum\limits_{k = 1}^n (X_k - \overline{X}_n)^2.

А теперь представим, что мы посчитали эти (или какие-то другие) статистики, а семейство распределений пока не выбрали. Как же совершить этот судьбоносный выбор? Давайте посмотрим на следующие три семейства и подумаем, в каком из них мы бы стали искать распределение, зная его истинные матожидание и дисперсию?

Three

Почему-то хочется сказать, что в первом. Почему? Второе не симметрично – но что нас может заставить подозревать, что интересующее нас распределение не симметрично? С третьим проблема в том, что, выбирая его, мы добавляем дополнительную информацию как минимум о том, что у распределения конечный носитель. А с чего бы? У нас такой инфомации вроде бы нет.

Общая идея такова: мы будем искать распределение, которое удовлетворяет только явно заданным нами ограничениям и не отражает никакого дополнительного знания о нём. Таким образом, искомое распределение должно обладать максимальной неопределённостью при заданных ограничениях, или, говоря более научно, иметь максимально возможную энтропию. В самом деле, энтропия выражает нашу меру незнания о том, как ведёт себя распределение, и чем она больше – тем более «‎произвольное» распределение, по крайней мере, в теории. В этом и заключается принцип максимальной энтропии для выбора модели машинного обучения.

Как мы уже видели выше, среди распределений с конечным носителем максимальную энтропию имеет равномерное распределение. Примеры геометрического и нормального распределения показывают, что энтропия распределений с бесконечным носителем (счётным или континуальным) может быть сколь угодно большой, и среди них нет какого-то одного распределения с максимальной энтропией. Однако в более узком классе распределений с фиксированным средним и/или дисперсией найти распределение с максимальной энтропией, как правило, можно.

Пример. Покажем, что среди распределений на множестве натуральных чисел N\mathbb N и математическим ожиданием μ>1\mu > 1 максимальную энтропию имеет геометрическое распределение.

Для минимизации энтропии H[p]=n=1pnlogpn\mathbb H[p] = -\sum\limits_{n=1}^\infty p_n \log p_n с учётом ограничений

n=1pn=1,n=1npn=μ\sum\limits_{n=1}^\infty p_n = 1, \quad \sum\limits_{n=1}^\infty np_n = \mu

воспользуемся методом множителей Лагранжа, согласно которому требуется минимизировать функцию Лагранжа

L(p,a,b)=n=1pnlogpna(n=1pn1)b(n=1npnμ). \mathcal L(p, a, b) = \sum\limits_{n=1}^\infty p_n \log p_n - a\Big(\sum\limits_{n=1}^\infty p_n - 1\Big) - b\Big(\sum\limits_{n=1}^\infty np_n - \mu\Big).

Приравняем к нулю частные производные по pnp_n:

Lpn=1+logpnabn=0. \frac{\partial \mathcal L}{\partial p_n} = 1 + \log p_n - a - bn = 0.

Отсюда следует, что pn=ea1+bn=αβnp_n = e^{a-1 + bn} = \alpha \beta^n, так что распределение действительно получается геометрическое. Параметры α\alpha и β\beta найдём из уравнений

1=n=1pn=n=1αβn=αβ1β, 1 = \sum\limits_{n=1}^\infty p_n = \sum\limits_{n=1}^\infty \alpha \beta^n = \frac{\alpha\beta}{1-\beta},

μ=n=1npn=αβn=1nβn1=αβ(1β)2. \mu = \sum\limits_{n=1}^\infty np_n = \alpha \beta\sum\limits_{n=1}^\infty n \beta^{n-1} = \frac{\alpha\beta}{(1-\beta)^2}.

Деля первое уравнение на второе, получаем 1μ=1β\frac 1 \mu = 1 - \beta, или β=11μ\beta = 1 - \frac 1\mu. Далее из первого уравнения находим α=1ββ=1μ1\alpha = \frac {1-\beta}{\beta} = \frac 1{\mu -1}. Итак,

pn=1μ1(11μ)n=1μ(11μ)n1,p_n = \frac 1{\mu -1} \Big(1 - \frac 1\mu\Big)^n = \frac 1\mu \Big(1 - \frac 1\mu\Big)^{n-1},

а это и есть геометрическое распределение с параметром 1μ\frac 1\mu.

У непрерывных распределений возможны более интересные комбинации из ограничений на носитель и параметры. И конечно же, первую скрипку среди распределений с максимальной энтропией играет гауссовское распределение.

Пример. Докажем, что среди распределений на R\mathbb R c математическим ожиданием μ\mu и дисперсией σ2\sigma^2 наибольшую энтропию имеет нормальное распределение N(μ,σ2)\mathcal{N}(\mu,\sigma^2).

Пусть p(x)p(x) – некоторое распределение со средним μ\mu и дисперсией σ2\sigma^2, q(x)N(μ,σ2)q(x)\sim\mathcal{N}(\mu, \sigma^2). Как было показано выше, H[q]=12log(2πσ2)+12\mathbb H[q] = \frac12\log(2\pi\sigma^2) + \frac12. Запишем дивергенцию Кульбака-Лейблера:

KL(pq)=+p(x)logp(x)dx+p(x)logq(x)dx=\mathbb KL(p\vert\vert q) = \int\limits_{-\infty}^{+\infty} p(x)\log{p(x)}dx - \int\limits_{-\infty}^{+\infty} p(x)\log{q(x)}dx =

=H[p]+p(x)(12log(2πσ2)12σ2(xμ)2)dx== -\mathbb H[p] - \int\limits_{-\infty}^{+\infty} p(x)\left(-\frac12\log(2\pi\sigma^2) - \frac 1{2\sigma^2}(x - \mu)^2\right)dx =

=H[p]+12log(2πσ2)+p(x)dx+12σ2+(xμ)2p(x)dx=V[p]=σ2== - \mathbb H[p] +\frac12\log(2\pi\sigma^2)\int\limits_{-\infty}^{+\infty} p(x)\,dx + \frac1{2\sigma^2}\underbrace{\int\limits_{-\infty}^{+\infty}(x - \mu)^2p(x)\,dx}_{=\mathbb{V}[p]=\sigma^2} =

=H[p]+12log(2πσ2)+12=H[q]H[p].= - \mathbb H[p] + \frac12\log(2\pi\sigma^2) + \frac12 = \mathbb H[q] - \mathbb H[p].

Так как KL-дивергенция всегда неотрицательна, получаем, что H[p]H[q]\mathbb H[p]\leqslant \mathbb H[q] при любом распределении pp, удовлетворяющем заданным ограничениям.

Можно показать, что максимальную энтропию среди многомерных распределений с вектором средних μ\boldsymbol \mu и матрицей ковариаций Σ\boldsymbol \Sigma имеет также гауссовское распределение N(μ,Σ)\mathcal N(\boldsymbol \mu, \boldsymbol \Sigma).

Упражнение. Докажите, что среди распределений на отрезке [a,b][a,b] максимальную энтропию имеет равномерное распределение U[a,b]U[a,b].

Попробуйте доказать самостоятельно, прежде чем смотреть решение.

Пусть q(x)U[a,b]q(x) \sim U[a,b], тогда для произвольного распределения pp на отрезке [a,b][a,b] имеем

H[p,q]=abp(x)logq(x)dx=abp(x)log(ba)dx=log(ba). \mathbb H[p,q] = - \int\limits_a^b p(x)\log q(x)\,dx = \int\limits_a^b p(x) \log(b-a)\,dx = \log(b-a).

В частности, отсюда следует, что H[q]=log(ba)\mathbb H[q] = \log(b-a). Расписывая теперь KL-дивергенцию, получаем

KL(pq)=H[p,q]H[p]=log(ba)H[p]=H[q]H[p]0,\mathbb{KL}(p\vert\vert q) = \mathbb H[p,q] - \mathbb H[p] = \log(b-a) - \mathbb H[p] = \mathbb H[q] - \mathbb H[p] \geqslant 0,

что и требовалось доказать.

Упражнение. Докажите, что среди распределений на промежутке [0,+)[0, +\infty) с математическим ожиданием λ>0\lambda > 0 максимальную энтропию имеет показательное распределение Exp(1λ)\mathrm{Exp\big(\frac 1\lambda)}.

Попробуйте доказать самостоятельно, прежде чем смотреть решение.

Согласно результату одного из предыдущих упражнений H[q]=logλ+1\mathbb H[q] = \log \lambda + 1, если q(x)Exp(1λ)q(x) \sim \mathrm{Exp}\big(\frac 1\lambda), т.е. q(x)=1λexλq(x) = \frac 1\lambda e^{-\frac x\lambda}, x0x\geqslant 0. Далее, для произвольного распределения p(x)p(x) на [0,+)[0, +\infty) со средним λ\lambda имеем

KL(pq)=H[p,q]H[p]=0+p(x)(logλxλ)dxH[p]=\mathbb{KL}(p\vert\vert q) = \mathbb H[p,q] - \mathbb H[p] = - \int\limits_0^{+\infty} p(x)\Big(-\log\lambda - \frac x\lambda\Big) \,dx - \mathbb H[p]=

=logλ+1λ0+xp(x)dxH[p]=logλ+1H[p]=H[q]H[p]0.=\log\lambda + \frac 1\lambda \int\limits_0^{+\infty} xp(x)\, dx - \mathbb H[p] = \log\lambda + 1 - \mathbb H[p] = \mathbb H[q]- \mathbb H[p] \geqslant 0.

Как выяснилось, многие классические распределения имеют максимальную энтропию при весьма естественных ограничениях. Но как быть, если даны не эти конкретные, а какие-то другие ограничения? Есть ли какой-нибудь надёжный алгоритм вывода распределения с максимальной энтропией, позволяющий избежать случайных озарений и гаданий на кофейной гуще? Оказывается, что при некоторых не очень обременительных ограничениях ответ можно записать с помощью распределений экспоненциального класса.

Экспоненциальное семейство распределений

Говорят, что параметрическое семейство распределений относится к экспоненциальному классу, если его pdf (или pmf) может быть представлена в виде

p(xθ)=g(x)h(θ)exp(θTu(x))=g(x)exp(θTu(x)A(θ)),p(x\vert \boldsymbol\theta) = \frac{g(x)}{h(\boldsymbol\theta)}\exp(\boldsymbol \theta^T \boldsymbol u(x)) = g(x)\exp(\boldsymbol \theta^T \boldsymbol u(x) - A(\boldsymbol \theta)),

где

  • θRm\boldsymbol\theta \in\mathbb R^m – вектор натуральных параметров распределения;
  • g(x)g(x) — неотрицательная функция (base measure), часто равная единице;
  • h(θ)>0h(\boldsymbol\theta) > 0нормализатор (partition), обеспечивающий суммируемость pmf или интегрируемость pdf в единицу:

h(θ)=g(x)exp(θTu(x))dx;h(\boldsymbol\theta) = \int g(x)\exp(\boldsymbol \theta^T \boldsymbol u(x))dx;

  • A(θ)=lnh(θ)A(\boldsymbol \theta) = \ln h(\boldsymbol\theta)log-partition;
  • u(x)Rm\boldsymbol u(x)\in\mathbb R^m — вектор достаточных статистик распределения.

Пример. Покажем, что нормальное распределение N(xμ,σ2)\mathcal N(x\vert \mu, \sigma^2) принадлежит экспоненциальному классу. Оно имеет два параметра, поэтому такую же размерность имеют θ\boldsymbol \theta и вектор-функция u\boldsymbol u.

Распишем плотность:

12πσexp((xμ)22σ2)=12πσexp(12σ2x2+μσ2x+μ22σ2)=exp(12σ2x2+μσ2x)2πσexp(μ22σ2).\frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right) = \frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac 1{2\sigma^2}x^2 + \frac{\mu}{\sigma^2}x + \frac{\mu^2}{2\sigma^2}\right)= \frac{\exp\left(-\frac1{2\sigma^2}x^2 + \frac{\mu}{\sigma^2}x\right)}{\sqrt{2\pi}\sigma\exp\left(-\frac{\mu^2}{2\sigma^2}\right)}.

Положим g(x)=12πg(x) = \frac 1{\sqrt{2\pi}}, u(x)=(x,x2)\boldsymbol u(x) = (x, x^2),

θ=(θ1,θ2),θ1=μσ2,θ2=12σ2.\boldsymbol \theta = (\theta_1, \theta_2),\quad \theta_1 = \frac{\mu}{\sigma^2},\quad \theta_2 = -\frac1{2\sigma^2}.

Остаётся выразить функцию h(θ)=σexp(μ22σ2)h(\boldsymbol\theta) = \sigma\exp\left(-\frac{\mu^2}{2\sigma^2}\right) через θ1\theta_1 и θ2\theta_2.

Упражнение. Выразите partition h(θ)h(\boldsymbol\theta) и log-partition A(θ)A(\boldsymbol\theta) через θ\boldsymbol\theta и запишите плотность нормального распределения в экспоненциальном виде.

Ответ

h(θ)=2θ2exp(θ124θ2),A(θ)=θ124θ212ln(2θ2),h(\boldsymbol\theta)=\sqrt{-2\theta_2} \exp\Big(\frac{\theta_1^2}{4\theta_2}\Big), \quad A(\boldsymbol\theta) = -\frac{\theta_1^2}{4\theta_2} - \frac 12\ln(-2\theta_2),

N(xμ,σ2)=12π1h(θ)exp(θTu(x))=12πexp(θTu(x)A(θ)).\mathcal N(x \vert \mu, \sigma^2) = \frac1{\sqrt {2\pi}}\cdot \frac 1{h(\boldsymbol\theta)}\exp(\boldsymbol\theta^T \boldsymbol u(x)) = \frac1{\sqrt {2\pi}} \exp(\boldsymbol \theta^T \boldsymbol u(x) - A(\boldsymbol \theta)).

Пример. Покажем, что распределение Бернулли Bern(p)\mathrm{Bern}(p) принадлежит экспоненциальному классу. Его pmf P(ξ=xp)\mathbb P(\xi =x \vert p) можно записать как

px(1p)1x=exp(xlogp+(1x)log(1p))=(1p)exp(xlogp1p). p^x(1 - p)^{1 - x} = \exp\big(x\log{p} + (1 - x)\log(1 - p)\big)= (1-p) \exp\Big(x\log\frac p{1-p}\Big).

Параметр здесь один, поэтому натуральный параметр θ\theta тоже один: θ=logp1p\theta = \log\frac p{1-p}. Такая функция от pp называется функцией логитов и активно участвует в построении модели логистической регрессии. Остальные функции положим равными u(x)=xu(x) = x, g(x)=1g(x)=1, h(θ)=11ph(\theta) = \frac 1{1-p}. Остаётся выразить partition через θ\theta:

logp1p=log(1+11p)=θ    11p=1+eθ. \log\frac p{1-p} = \log\Big(-1 + \frac 1{1-p}\Big) = \theta \iff \frac 1{1-p} = 1 + e^\theta.

Итак, h(θ)=1+eθh(\theta) = 1+e^\theta, и экспоненциальный вид распределения Бернулли записывается как

11+eθexp(θx)=exp(θxlog(1+eθ)). \frac 1{1+e^\theta}\exp(\theta x) = \exp\big(\theta x - \log(1 + e^\theta)\big).

Вопрос на подумать. Принадлежит ли к экспоненциальному классу семейство равномерных распределений на отрезках U[a,b]U[a, b]? Казалось бы, да, ведь

p(x)=1baI[a,b](x)exp(0).p(x) = \frac{1}{b - a}\mathbb{I}_{[a,b]}(x)\exp(0).

В чём может быть подвох?

Попробуйте определить сами, прежде чем смотреть ответ.

Нет, не принадлежит. Давайте вспомним, как звучало определение экспоненциального семейства. Возможно, вас удивило, что там было написано не «распределение относится», а «семейство распределений относится». Это важно: ведь семейство определяется именно различными значениями θ\theta, и если нас интересует семейство равномерных распределений на отрезках, определяемое параметрами aa и bb, то они не могут быть в функции g(x)g(x), они должны быть под экспонентой, а экспонента ни от чего не может быть равна индикатору.

При этом странное и не очень полезное семейство с нулём параметров, состоящее из одинокого распределения U[0,1]U[0,1] можно считать относящимся к экспоненциальному классу: ведь для него формула

p(x)=I[0,1](x)exp(0)p(x) = \mathbb{I}_{[0,1]}(x)\exp(0)

будет работать.

К экспоненциальным семействам относятся многие непрерывные и дискретные распределения из часто встречающихся в теории и на практике, в том числе

  • нормальное N(μ,σ2)\mathcal N(\mu, \sigma^2);
  • распределение Пуассона Pois(λ)\mathrm{Pois}(\lambda);
  • экспоненциальное Exp(λ)\mathrm{Exp}(\lambda);
  • биномиальное Bin(n,p)\mathrm{Bin}(n, p);
  • геометрическое Geom(p)\mathrm{Geom}(p);
  • бета-распределение;
  • гамма-распределение;
  • распределение Дирихле.

Как выглядят натуральные параметры, достаточные статистики и нормализаторы этих и других распределений из экспоненциального класса, можно посмотреть на википедии.

К экспоненциальным семействам не относятся, например, равномерное распределение U[a,b]U[a,b], tt-распределение Стьюдента, распределение Коши, смесь нормальных распределений.

Дифференцирование log-partition

Если распределение p(xθ)p(x\vert \boldsymbol\theta) принадлежит экспоненциальному классу,

p(xθ)=g(x)h(θ)exp(θTu(x))=g(x)exp(θTu(x)A(θ)),p(x\vert \boldsymbol\theta) = \frac{g(x)}{h(\boldsymbol\theta)}\exp(\boldsymbol \theta^T \boldsymbol u(x)) = g(x)\exp(\boldsymbol \theta^T \boldsymbol u(x) - A(\boldsymbol \theta)),

то моменты его достаточных статистик u(x)\boldsymbol u(x) могут быть получены дифференцированием функции A(θ)=logh(θ)A(\boldsymbol \theta) = \log h(\boldsymbol\theta).

Утверждение. θlogh(θ)=Eξp(xθ)u(ξ)\nabla_{\boldsymbol\theta}\log h(\boldsymbol \theta) = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi).

Доказательство.
По правилу дифференцирования сложной функции имеем

θlogh(θ)=θh(θ)h(θ). \nabla_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} = \frac{\nabla_{\boldsymbol\theta}{h(\boldsymbol\theta)}}{h(\boldsymbol\theta)}.

Нормализатор h(θ)h(\boldsymbol\theta) записывается в виде интеграла

h(θ)=g(x)exp(θTu(x))dx, h(\boldsymbol\theta) = \int g(x)\exp\left(\boldsymbol\theta^T\boldsymbol u(x)\right)dx,

который мы продифференцируем внесением градиента внутрь под знак интеграла:

θh(θ)=θg(x)exp(θTu(x))dx=\nabla_{\boldsymbol\theta} h(\boldsymbol\theta) = \nabla_{\boldsymbol\theta} \int g(x)\exp\left(\boldsymbol\theta^T\boldsymbol u(x)\right)dx =

=g(x)θexp(θTu(x))dx=g(x)u(x)exp(θTu(x))dx. =\int g(x)\nabla_{\boldsymbol\theta} \exp\left(\boldsymbol\theta^T\boldsymbol u(x)\right)dx = \int g(x)\boldsymbol u(x) \exp\big(\boldsymbol\theta^T\boldsymbol u(x)\big)dx.

Таким образом,

θh(θ)h(θ)=u(x)g(x)h(θ)exp(θTu(x))p(xθ)dx=Eξp(xθ)u(ξ). \frac{\nabla_{\boldsymbol\theta}{h(\boldsymbol\theta)}}{h(\boldsymbol\theta)} =\int \boldsymbol u(x)\underbrace{\frac{g(x)}{h(\boldsymbol \theta)}\exp(\boldsymbol \theta^T\boldsymbol u(x))}_{p(x\vert\boldsymbol \theta)}dx = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi).

Замечание о дифференцировании под знаком интегралом

Использованный в доказательстве приём внесения градиента под знак интеграла называют также правилом Лейбница. Этот же метод используется для почленного дифференцирования ряда, что может быть полезно в случае дискретного распределения p(xθ)p(x\vert \boldsymbol\theta). В математическом анализе имеется ряд теорем, обеспечивающих применимость правила Лейбница, однако, мы не будем на них останавливаться. Будем считать, что все рассматриваемые экспоненциальные семейства таковы, что применение правила Лейбница законно.

Если ui(x)=xiu_i(x) = x^i, то в соответствии с только что доказанным частная производная A(θ)θi\frac{\partial A(\boldsymbol \theta)}{\partial \theta_i} даёт ii-й момент распределения p(xθ)p(x\vert \boldsymbol\theta).

Упражнение. Вычислите производные по натуральным параметрам от log-partition для распределения Бернулли Bern(xp)\mathrm{Bern}(x\vert p) и нормального распределения N(μ,σ2)\mathcal N(\mu, \sigma^2) и проверьте, что они совпадают со значениями соответствующих моментов.

Решение

Для ξBern(xp)\xi \sim \mathrm{Bern}(x\vert p) имеем

A(θ)=log(1+eθ),11p=1+eθ, A(\theta) = \log(1 + e^\theta), \quad \frac 1{1-p} = 1 + e^\theta,

поэтому

A(θ)=eθ1+eθ=111+eθ=11+p=p=Eξ. A'(\theta) = \frac{e^\theta}{1+e^\theta} = 1 - \frac 1{1+e^\theta} = 1 -1 + p = p = \mathbb E\xi.

Для гауссовской случайной величины ηN(xμ,σ2)\eta\sim\mathcal N(x \vert\mu, \sigma^2) получаем

A(θ1,θ2)=θ124θ212ln(2θ2), A(\theta_1, \theta_2) = -\frac{\theta_1^2}{4\theta_2} - \frac 12\ln(-2\theta_2),

θ1=μσ2, \quad \theta_1 = \frac{\mu}{\sigma^2},

θ2=12σ2, \quad \theta_2 = -\frac1{2\sigma^2},

и, значит,

Aθ1=θ12θ2=μ=Eη, \frac{\partial A}{\partial\theta_1} = -\frac{\theta_1}{2\theta_2} = \mu = \mathbb E\eta,

Aθ2=θ124θ2212θ2=μ2+σ2=Eη2. \frac{\partial A}{\partial\theta_2} = \frac{\theta_1^2}{4\theta_2^2} - \frac 1{2\theta_2} = \mu^2 + \sigma^2 = \mathbb E\eta^2.

Кстати, можно продифференцировать ещё раз и доказать, что

θ2logh(θ)=cov(u(ξ),u(ξ)). \nabla^2_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} = \mathrm{cov}(\boldsymbol u(\xi), \boldsymbol u(\xi)).

Доказательство

В предыдущий раз мы доказали, что

θA(θ)=u(x)g(x)exp(θTu(x)A(θ))dx.\nabla_{\boldsymbol\theta}A(\boldsymbol\theta) = \int \boldsymbol u(x)g(x)\exp(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta))\,dx.

Теперь возьмём ещё раз градиент по θ\boldsymbol \theta от этого выражения:

θ2A(θ)=u(x)g(x)θexp(θTu(x)A(θ))dx=\nabla^2_{\boldsymbol\theta}A(\boldsymbol\theta) = \int \boldsymbol u(x)g(x)\nabla_{\boldsymbol\theta}\exp\big(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)\big)\,dx =

=u(x)g(x)(u(x)θA(θ))Texp(θTu(x)A(θ))dx== \int \boldsymbol u(x)g(x)\big(\boldsymbol u(x) - \nabla_{\boldsymbol\theta}A(\boldsymbol \theta)\big)^T \exp\big(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)\big)\,dx=

=u(x)g(x)(u(x)Eu)Texp(θTu(x)A(θ))dx==\int \boldsymbol u(x)g(x)\big(\boldsymbol u(x) - \mathbb{E} \boldsymbol u\big)^T \exp\big(\boldsymbol \theta^T\boldsymbol u(x) - A(\boldsymbol \theta)\big)\,dx =

=EuuTEu(Eu)T=cov(u(ξ),u(ξ)).= \mathbb{E} \boldsymbol u \boldsymbol u^T - \mathbb{E} \boldsymbol u \big(\mathbb{E} \boldsymbol u\big)^T = \mathrm{cov}(\boldsymbol u(\xi), \boldsymbol u(\xi)).

В последней выкладке для краткости обозначено Eu=Eξp(xθ)u(ξ)\mathbb E \boldsymbol u = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi).

MLE для семейства из экспоненциального класса

Возможно, вас удивил странный и на первый взгляд не очень естественный вид p(xθ)p(x\vert \boldsymbol \theta). Но всё не просто так: оказывается, что оценка максимального правдоподобия параметров распределений из экспоненциального класса устроена очень интригующе.

Запишем функцию правдоподобия i.i.d. выборки x1,,xnx_1,\ldots,x_n:

p(x1,,xnθ)=h(θ)n(i=1ng(xi))exp(θTi=1nu(xi)).p(x_1, \ldots, x_n\vert\boldsymbol \theta) = h(\boldsymbol \theta)^{-n}\cdot\bigg(\prod_{i=1}^ng(x_i)\bigg)\cdot\exp\bigg(\boldsymbol \theta^T \sum_{i=1}^n \boldsymbol u(x_i)\bigg).

Её логарифм равен

logp(x1,,xnθ)=nlogh(θ)+i=1nlogg(xi)+θTi=1nu(xi).\log p(x_1, \ldots, x_n\vert\boldsymbol \theta) = -n\log{h(\boldsymbol \theta)} + \sum_{i=1}^n\log{g(x_i)} + \boldsymbol \theta^T\sum_{i=1}^n \boldsymbol u(x_i).

Дифференцируя по θ\boldsymbol \theta, получаем

θlogp(x1,,xnθ)=nθlogh(θ)+i=1nu(xi).\nabla_{\boldsymbol\theta}\log p(x_1, \ldots, x_n\vert\boldsymbol\theta) = -n\nabla_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} + \sum_{i=1}^n \boldsymbol u(x_i).

Приравнивая θlogp(x1,,xnθ)\nabla_{\boldsymbol\theta}\log p(x_1,\ldots, x_n\vert\boldsymbol\theta) к нулю и пользуясь равенством θlogh(θ)=Eξp(xθ)u(ξ)\nabla_{\boldsymbol\theta}\log{h(\boldsymbol\theta)} = \mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi), находим

Eξp(xθ)u(ξ)=1ni=1nu(xi).\mathbb{E}_{\xi \sim p(x\vert \boldsymbol \theta)} \boldsymbol u(\xi) = \frac1n \sum_{i=1}^n \boldsymbol u(x_i).

Таким образом, теоретические матожидания всех компонент ui(ξ)u_i(\xi) должны совпадать с их эмпирическими оценками, а метод максимального правдоподобия совпадает с методом моментов для Eui(ξ)\mathbb{E}u_i(\xi) в качестве моментов. И в следующем пункте выяснится, что распределения из экспоненциальных семейств обладают максимальной энтропией среди тех, что имеют заданные моменты Eui(ξ)\mathbb{E}u_i(\xi).

Теорема Купмана—Питмана—Дармуа

Теперь мы наконец готовы сформулировать одно из самых любопытных свойств семейств экспоненциального класса.

В следующей теореме мы опустим некоторые не очень обременительные условия регулярности. Просто считайте, что для хороших дискретных и абсолютно непрерывных распределений, с которыми вы в основном и будете сталкиваться, это так.

Теорема. Пусть параметр θRm\boldsymbol\theta\in R^m распределения pθ(x)=1h(θ)exp(θTu(x))p_{\boldsymbol\theta}(x) = \frac1{h(\boldsymbol\theta)}\exp\big(\boldsymbol\theta^T\boldsymbol u(x)\big) выбран так, что

Eξpθ(x)u(ξ)=α\mathbb{E}_{\xi \sim p_{\boldsymbol\theta}(x)} \boldsymbol u(\xi) = \boldsymbol\alpha

для некоторого фиксированного αRm\boldsymbol \alpha \in\mathbb R^m. Тогда распределение pθ(x)p_{\boldsymbol\theta}(x) обладает наибольшей энтропией среди распределений qq с тем же носителем, для которых Eξq(x)u(ξ)=α\mathbb{E}_{\xi \sim q(x)} \boldsymbol u(\xi) = \boldsymbol\alpha.

Идея обоснования через оптимизацию.

Мы приведём рассуждение для дискретного случая; в абсолютно непрерывном рассуждения будут по сути теми же, только там придётся дифференцировать не по переменных, а по функциям, и мы решили не ввергать читателя в мир вариационного исчисления.

В дискретном случае у нас есть счётное семейство точек x1,x2,x_1, x_2,\ldots, и распределение определяется счётным набором вероятностей pip_i принимать значение xix_i. Мы будем решать задачу

{jpjlogpjmax,jpjui(xj)=αi,i=1,,m,jpj=1,pj0.\begin{cases} -\sum_j p_j\log{p_j}\longrightarrow\max,\\ \sum_jp_ju_i(x_j) = \alpha_i, i = 1,\ldots,m,\\ \sum_jp_j = 1,\\ p_j\geqslant 0. \end{cases}

Для решения этой оптимизационной задачи нам понадобится обобщение метода множителей Лагранжа, известное также как теорема Каруша—Куна—Таккера. В данном случае задача сводится к минимизации лагранжиана

L=jpjlogpj+iθi(αijpjui(xj))+θ0(jpj1)jλjpj\mathcal{L} = \sum_j p_j\log{p_j} + \sum_i\theta_i\left(\alpha_i - \sum_jp_ju_i(x_j)\right)+\theta_0\left(\sum_jp_j - 1\right) - \sum_j\lambda_jp_j

при имеющихся ограничениях и условиях дополняющей нежёсткости λjpj=0\lambda_jp_j = 0.
Приравняем частные производные по pjp_j к нулю:

Lpj=logpj+1iθiui(xj)+θ0λj=0.\frac{\partial\mathcal{L}}{\partial p_j} = \log{p_j} + 1 - \sum_i\theta_iu_i(x_j) + \theta_0 - \lambda_j = 0.

Отсюда получаем, что

pj=exp(θTu(xj))exp(λjθ01).p_j = \frac{\exp(\boldsymbol \theta^T \boldsymbol u(x_j))}{\exp\left(\lambda_j - \theta_0 - 1\right)}.

Числитель уже ровно такой, как и должен быть у распределения из экспоненциального класса; разберёмся со знаменателем. Поскольку pj>0p_j >0 (ведь тут сплошные экспоненты), λj=0\lambda_j = 0 при всех jj. Параметр θ0\theta_0 находится из условия jpj=1\sum_jp_j = 1, а точнее, выражается через остальные θi\theta_i, что позволяет записать знаменатель в виде h(θ)h(\boldsymbol \theta).

Идея доказательства «в лоб».

Как и следовало ожидать, оно ничем не отличается от того, как мы доказывали максимальность энтропии у равномерного или нормального распределения. Пусть q(x)q(x) – ещё одно распределение, для которого

ui(x)q(x)dx=ui(x)pθ(x)dx\int u_i(x)q(x)dx = \int u_i(x)p_{\boldsymbol\theta}(x)dx

для всех i=1,,mi = 1,\ldots,m. Тогда

0KL(qp)=q(x)log(q(x)pθ(x))dx=0\leqslant \mathbb{KL}(q\vert\vert p) = \int q(x)\log\left(\frac{q(x)}{p_{\boldsymbol\theta}(x)}\right)dx =

=q(x)logq(x)dxH[q]q(x)logpθ(x)dx==\underbrace{\int q(x)\log{q(x)}dx}_{-\mathbb H[q]} - \int q(x)\log{p_{\boldsymbol\theta}(x)}dx=

=H[q]q(x)(logh(θ)+iθiui(x))dx==-\mathbb H[q] - \int q(x)\left(-\log{h(\boldsymbol \theta)} + \sum_i\theta_iu_i(x)\right)dx =

=H[q]logh(θ)q(x)dx=1=pθ(x)dxiθiq(x)ui(x)dx=pθ(x)ui(x)dx==-\mathbb H[q] - \log{h(\theta)}\underbrace{\int q(x)dx}_{=1=\int p_{\boldsymbol\theta}(x)dx} - \sum_i\theta_i\underbrace{\int q(x)u_i(x)dx}_{=\int p_{\boldsymbol\theta}(x)u_i(x)dx} =

=H[q]p(x)(logh(θ)+iθiui(x))dx==-\mathbb H[q] - \int p(x)\left(-\log{h(\boldsymbol\theta)} + \sum_i\theta_iu_i(x)\right)dx =

=H[q]+pθ(x)logpθ(x)dx=H[q]+H[pθ].=-\mathbb H[q] + \int p_{\boldsymbol\theta}(x)\log{p_{\boldsymbol\theta}(x)}dx = -\mathbb H[q] + \mathbb H[p_{\boldsymbol\theta}].

Таким образом, H[pθ]H[q]\mathbb H[p_{\boldsymbol\theta}]\geqslant \mathbb H[q], что и требовалось доказать.

Выше мы уже находили обладающее максимальной энтропией распределение на множестве натуральных чисел с заданным математическим ожиданием μ>1\mu>1. Таковым оказалось геометрическое распределение Geom(1μ)\mathrm{Geom}\big(\frac 1\mu\big).

Теорема Купмана—Питмана—Дармуа позволяет сделать это гораздо быстрее.
В данном случае у нас лишь одна функция u1(x)=xu_1(x) = x, которая соответствует фиксации математического ожидания Eξ\mathbb E\xi. Искомое дискретное распределение имеет вид

pk=P(ξ=k)=1h(θ)exp(θk)=1h(θ)(eθ)k.p_k =\mathbb P(\xi = k) = \frac1{h(\theta)}\exp(\theta k) = \frac1{h(\theta)} \big(e^{\theta}\big)^k.

Это уже похоже на геометрическое распределение с параметром p=1eθp = 1 - e^{\theta}. Его математическое ожидание равно 1p\frac 1p, что по условию должно равняться μ\mu. Итак, наше распределение с максимальной этропией выглядит так:

pk=1μ(11μ)k1,kN.p_k = \frac1{\mu}\left(1 - \frac1{\mu}\right)^{k-1},\quad k\in\mathbb N.

Пример. Среди распределений на всей вещественной прямой с заданным математическим ожиданием μ\mu найдём распределение с максимальной энтропией.

А сможете ли вы его найти? Решение под катом.

Теория говорит нам, что его плотность должна иметь вид

p(x)=1h(θ)exp(θx),p(x) = \frac1{h(\theta)}\exp\left(\theta x\right),

но интеграла экспоненты не существует, то есть применение «в лоб» теоремы провалилось. И неспроста: если даже рассмотреть все нормально распределённые случайные величины со средним μ\mu, их энтропии, равные 12+12log(2πσ2)\frac12 + \frac12\log(2\pi\sigma^2), не ограничены сверху, то есть величины с наибольшей энтропией не существует.

Чтобы добавить в заметки выделенный текст, нажмите Command + E

Пройдите квиз по параграфу

Чтобы закрепить пройденный материал
Предыдущий параграф16.5. Независимость и условные распределения вероятностей
Предыдущий параграф16.6. Параметрические оценки

Вступайте в сообщество хендбука

Здесь можно найти единомышленников, экспертов и просто интересных собеседников. А ещё — получить помощь или поделиться знаниями.